%-------------------------------------------------------------------------------
% $Id: wav1d_sep_vars.m,v 1.3 2011/06/20 21:53:39 paul Exp $
% $Date: 2011/06/20 21:53:39 $
% $Author: paul $
%
%% Separation of Variables Solution to 1D wave equation of a vibrating string
%-------------------------------------------------------------------------------
N = 1024;                      % number of samples on x-axis
M = 500;                       % number of samples in time interval T
T = 2;                         % time interval (0,T)

total_length = 1.0;            % length of string on x-axis
disturbance_length = 0.2;      % interval of disturbance funtion is (0,0.2)
                               
dx = 2 * total_length/N;
dt = T / M;

N1 = floor(0.2/dx);
N2 = floor((2*total_length-0.2)/dx + 2);

xs = [0:dx:2*total_length];
fs = zeros(1,N+1);
fs(2:N1) = exp(-1./(1-(10*xs(2:N1)-1).^2))
fs(N:-1:N2+1) = -fs(2:N1);

alpha = fft(fs(1:N)/N);
b = real(2 * i * alpha(1:N/2));

c = 1.0;                      

%-------------------------------------------------------------------------------
% Setup AVI file and parameters
aviobj = avifile('wav1d_sep_vars.avi');
aviobj.compression = 'None';             % no other option for Unix
aviobj.fps = 5;                          % frames per second
aviobj.quality = 10;                     % low quality is OK

%-------------------------------------------------------------------------------
% Animation loop
fig = figure;

for time=0:dt:T
  ss = zeros(1,N/2);
  for k=1:N/2-1
    ss = ss + b(k+1)*sin(k*pi*xs(1:N/2)/total_length) * cos(k*pi*time/total_length);
  end
  
  plot([1:N/2],ss);
  axis([0 512 -0.4 0.4]);
  title('1D Wave Equation Separation of Variables Solution');
  frame = getframe(fig); 
  aviobj = addframe(aviobj,frame);
end 

%-------------------------------------------------------------------------------
% Clean up rendering frame and AVI file

close(fig);
aviobj = close(aviobj);